Anti-colorectal cancer of Ardisia gigantifolia Stapf. and targets prediction via network pharmacology and molecular docking study

Background Ardisia gigantifolia Stapf. (AGS), a Chinese folk medicine widely grows in the south of China and several studies reported that AGS could inhibit the proliferation of breast cancer, liver cancer, and bladder cancer cell lines. However, little is known about its anti-colorectal cancer (CRC) efficiency. Methods In the present study, a combination of MTT assay, network pharmacological analysis, bioinformatics, molecular docking, and molecular dynamics simulation study was used to investigate the active ingredients, and targets of AGS against CRC, as well as the potential mechanism. Results MTT assay showed that three kinds of fractions from AGS, including the n-butanol extract (NBAGS), ethyl acetate fraction (EAAGS), and petroleum ether fraction (PEAGS) significantly inhibited the proliferation of CRC cells, with the IC50 values of 197.24, 264.85, 15.45 µg/mL on HCT116 cells, and 523.6, 323.59, 150.31 µg/mL on SW620 cells, respectively. Eleven active ingredients, including, 11-O-galloylbergenin, 11-O-protocatechuoylbergenin, 11-O-syringylbergenin, ardisiacrispin B, bergenin, epicatechin-3-gallate, gallic acid, quercetin, stigmasterol, stigmasterol-3-o-β-D-glucopyranoside were identified. A total of 173 targets related to the bioactive components and 21,572 targets related to CRC were picked out through database searching. Based on the crossover targets of AGS and CRC, a protein-protein interaction network was built up by the String database, from which it was concluded that the core targets would be SRC, MAPK1, ESR1, HSP90AA1, MAPK8. Besides, GO analysis showed that the numbers of biological process, cellular component, and molecular function of AGS against CRC were 1079, 44, and 132, respectively, and KEGG pathway enrichment indicated that 96 signaling pathways in all would probably be involved in AGS against CRC, among which MAPK signaling pathway, lipid, and atherosclerosis, proteoglycans in cancer, prostate cancer, adherens junction would probably be the major pathways. The docking study verified that AGS had multiple ingredients and multiple targets against CRC. Molecular dynamics (MD) simulation analysis showed that the binding would be stable via forming hydrogen bonds. Conclusion Our study showed that AGS had good anti-CRC potency with the characteristics of multi-ingredients, -targets, and -signaling pathways. Supplementary Information The online version contains supplementary material available at 10.1186/s12906-022-03822-8.


Background
Colorectal cancer (CRC) is the third most commonly diagnosed cancer worldwide which occurs in the colon or rectum [1]. The incidence of CRC dramatically increased in the past few decades because of the changes in human lifestyle, environment, and aged populations [2]. In 2018, the International Agency for Research on Cancer (IARC) reported that about 1.8 million new cases of CRC were diagnosed (approximately 10.2% of total cancer cases) and 860,000 CRC-associated deaths occurred worldwide (approximately 9.2% of all cancer-related deaths) [3]. The new cases and CRC-associated deaths were predicted to increase to 2.2 million and 1.1 million by 2030, respectively [4]. The etiology of CRC is highly complicated and some CRC are genetically predisposed such as familial adenomatous polyposis (FAP), hereditary nonpolyposis colorectal cancer (HNPCC or Lynch syndrome), Peutz-Jeghers syndrome, and so on [5,6]. Some modifiable risk factors, including smoking, western-style diet addiction, obesity, diabetes, alcohol over-consumption, physical inactivity, antibiotic abuse, and intestinal microbiota disorder, were reported to play an important role in the pathogenesis of CRC [7,8].
In the clinic, the conventional treatments of CRC mainly include surgery, chemotherapy/radiotherapy, and targeted drug therapy. The CRC drugs currently used are 5-fluorouracil (5-FU), oxaliplatin, and irinotecan, and monoclonal antibodies newly developed including bevacizumab (a monoclonal antibody against VEGF) and cetuximab (a monoclonal antibody against EGFR). However, the 5-year survival rate of CRC patients in the advanced stages remains poor, being 18.5% in the United States and 27.7% in Europe, respectively [3]. Hence, more effective treatments or alternative remedies for CRC are urgently needed. Traditional Chinese herbal medicine (TCM) is a big resource for new drug development due to its long-time usage in the clinic, and numbers of cancer cases proved its good efficiency in neoplasia prevention. A multi-center prospective cohort study suggested that a long-time usage of TCM increased the survival rate of CRC patients at stages II and III [9]. A retrospective cohort study also indicated that TCM could significantly improve disease-free survival, in particular for patients with stage III CRC [10].
Ardisia gigantifolia Stapf. (AGS) is a kind of TCM that widely grows in the south of China. It was first recorded in the book of Sheng-Cao-Yao-Xing-Bei-Yao in the Qing dynasty and was commonly used for treating rheumatism, the pain of bones and muscles, and traumatic injury [11]. Previous studies showed that AGS mainly contained phenols, quinones, sterols, coumarins, triterpenoids, volatile oils, and flavonoids [12,13]. AGS had a wide range of biological properties such as antiinflammation, -oxidation, -thrombosis, and -cancer [14]. Several studies reported that AGS could inhibit the proliferation of MDA-MB-231 cells (breast cancer) [15], Bel-7402 cells (liver cancer) [16], and EJ cells (bladder cancer) [17]. However, there is no study reporting the anti-CRC efficiency of AGS and its underlying mechanism. This study is designed to investigate the anti-proliferation of AGS against CRC and explore the potential molecular mechanism.
Network pharmacology, an emerging methodology and a useful bioinformatics tool to investigate the complex effects and mechanisms between drugs and diseases at molecular, cellular, tissue, and biologic levels from a systems-level perspective, is now widely used in TCM investigation [18]. In this study, the growth inhibition of AGS against CRC cells was investigated via MTT assay. Network pharmacology was used to analyze the active ingredients, and potential targets of AGS in anti-CRC, as well as predict the possible mechanism and signaling pathway. Furthermore, molecular docking and molecular dynamics simulation were performed to study the binding pattern and stability between active compounds of AGS and therapeutic targets of CRC. The whole work procedure is shown in Fig. 1.

Samples preparation of AGS
The roots and rhizomes of AGS were powdered by a high-speed multifunctional crusher (JP-300 A-8, Yongkang, China) to obtain a granulometry of 2 ~ 5 mm. 3 kg of dry powder was combined with 70% ethanol in a round glass flask and extracted by a reflux method for 2 h. This process was repeated three times with the ratio of powder/ethanol of 1:8, 1:6, and 1:4 (w/v), respectively. All the extracted products were collected and dried under a vacuum freeze dryer (LGJ-18, Beijing, China). The extraction yield is 11.28%. Afterward, the ethanolic extracts were dissolved in distilled water and mixed well with diatomite. Dry it in an automatic program-controlled oven with the temperature setting at 70 ℃ and then grind the mixture into a fine powder [19]. The mixture was sequentially macerated and extracted using petroleum ether, ethyl acetate, and n-butanol, respectively, to obtain the three parts of extracts. The extracts were concentrated via a rotary evaporator and lyophilized to obtain three kinds of dry powder, naming PEAGS (, 4.25 mg/g of dried AGS), EAAGS (ethyl acetate extract, 5.1 mg/g of dried AGS), and NBAGS (n-butanol extract, 17.1 mg/g of dried AGS). The extracts were stored at − 35 °C until use.

Cells culture
Human colorectal cancer cell lines (HCT-116 and SW620) were commercially purchased from the cell resource center of the Shanghai Institutes for Biological Sciences, Chinese Academy of Sciences (Shanghai, China), and cultured in McCoy's 5 A medium and Leibovitz's L-15 medium supplemented with 10% fetal bovine serum and 1% penicillin-streptomycin, respectively, at 37 ℃ in a humidified incubator with 5% CO 2 .

MTT assay
HCT-116 cells (2.5 × 10 3 / well) and SW620 cells (4.5 × 10 3 / well) were seeded in 96-well plates and allowed to grow for 24 h, and then treated with different concentrations of PEAGS, EAAGS, NBAGS (0 ~ 1000 µg/ mL) for 48 h. After treatment, 20 µL of MTT solution (5 mg/mL) was added to each well and co-incubated for another 4 h at 37 °C. Finally, the supernatant was removed and 150 µL of DMSO was added to each well to dissolve the purple formazan crystals. A microplate reader (PerkinElmer, USA) was used to measure the optical density (OD) of each well at a wavelength of 570 nm. Each experiment was performed in triplicate.

Network pharmacology analysis
The potential active ingredients in AGS were determined by referring to the published literature by searching 3 databases, the CNKI database (https:// www. cnki. net/), Blyun database (http:// www. blyun. com/), and Pubmed database (https:// pubmed. ncbi. nlm. nih. gov/). The chemical structures of these active ingredients were drawn using the ChemDraw software and were input into the PharmMapper database (http:// www. lilab-ecust. cn/ pharm mapper/) to obtain the PDB IDs, target names, and fit scores. Then, the potential target set was further summarized according to the fit scores (> 0.7). After excluding the same targets and non-Homo sapiens targets from the screened targets, the PDB IDs were transformed into the gene symbol and gene IDs via the UniProt database (http:// www. unipr ot. org/). The CRC-related targets were obtained from the GeneCards Human database (https:// www. genec ards. org) using "colorectal cancer" as a keyword. The crossover genes between AGS and CRC were screened by the R software using the Venn Diagram package after eliminating duplicates. A compound-disease-target network was constructed using the Cytoscape3.7.2 to further explore the therapeutic mechanism of AGS against CRC. The PPI network construction of the predicted targets of AGS in treating CRC was performed based on their interaction data by applying the Interacting Genes/Proteins (STRING) database (version 11.0, https:// string-db. org/) and visualized using Cytoscape v3.7.2 with the parameters of the minimum required interaction score = 0.9 and hide disconnected nodes in the network. The network can be used to predict the protein interactions (including the physical and functional association between protein targets), among which the nodes represent the intersected target proteins, and edges represent the predicted or validated interactions between proteins. All collective proteins/genes were subjected to Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis using the database (DAVID, http:// david. abcc. ncifc rf. gov/) for annotation, visualization, and integrated discovery. GO analysis consists of Biological Process (BP), Cellular Component (CC), and Molecular Function (MF). KEGG analysis [20][21][22] aims to identify the significantly altered metabolic pathways of AGS against CRC based on a bioinformatics resource. GO analysis (p < 0.05) and KEGG enrichment analysis (p < 0.05) were visualized using the bioinformatics platform (http:// www. bioin forma tics. com. cn/).

Molecular docking analysis
To explore the interactions and binding modes between the active ingredients of AGS and the predicted targets from the network analysis, a molecular docking  [23], MAPK1 (PBD ID: 1PME) [24], ESR1 (PBD ID: 1A52) [25], HSP90AA1 (PBD ID: 7lt0) [26], MAPK8 (PBD ID: 1UKI) [27] were downloaded from the RCSB Protein Data Bank (http:// www. rcsb. org). All the proteins were processed by undergoing the following procedure before the docking calculations: preprocessing, reviewing, modification, and refinement (including optimization, removing waters, and minimizing the energy using the OPLS3e force field). The binding sites were determined according to the coordinates of the protein ligands. All compounds were docked by docking with standard precision (SP). The docked conformers were evaluated using the docking score. PyMOL 2.4.0 software was used to study the interaction between the docked molecules.

Molecular dynamic simulation
Molecular dynamics (MD) simulation analysis was performed to evaluate the stability of the protein-small molecule complex in the present study. MD can not only calculate the binding affinities of small molecules within the binding sites in the binding process but also exhibit dynamic conformational changes with the time scale [28,29]. Therefore, 100 ns of MD simulation was conducted and a total of 10 protein-small molecule complexes were subjected to MD simulations, including SRC-11-O-protocatechuoylbergenin complex, SRC-epicatechin-3-gallate complex, MAPK1-11-Ogalloylbergenin complex, MAPK1-epicatechin-3-gallate complex, ESR1-bergenin complex, ESR1-epicatechin-3-gallate complex, HSP90AA1-11-O-protocatechuoylbergenin complex, HSP90AA1-11-O-syringylbergenin complex, MAPK8-11-O-galloylbergenin complex. For the simulation, the AMBER18 package [30] was used to prepare and equilibrate the system with AMBER Force Field ff14SB for the proteins and Force Field GAFF2 for the small molecules [31,32]. Before the simulation, partial charges for small molecules were calculated by using the antechamber module and Gaussian at the Hartree-Fock (HF) SCF/6-31G* level of theory [33,34]. A rectangular periodic box of pre-equilibrated threepoint transferable intermolecular potential (TIP3P) solvent was used with a minimum distance of 10 Å [35] for solvating the complexes. The electroneutrality of protein/small molecule systems was maintained by adding an appropriate amount of sodium (NA) and chloride (CL) ions if needed. In the next step, the energy minimization (EM) was executed by 2500 steps of steepest descent followed by 2500 steps of a conjugate gradient. After energy minimization, each system was gradually heated from 0 to 298.15 K in a double time of 100 ps with position restraints. After that, the equilibration with position restraint on the protein was performed for 500 ps using NVT (number of particles, volume, and temperature) and NPT (number of particles, pressure, and temperature) ensembles with the temperature of 298.15 K and pressure of 1 bar, respectively. The Particle Mesh Ewald (PME) method was applied to calculate all the long-range electrostatic interactions during the MD simulations process with a radius of 10 Å for coulomb interactions [36]. SHAKE algorithm [37] and Langevin method [38] were performed to constrain all bonds and control the temperature, respectively. Finally, a 100 ns MD was performed with a time step of 2 fs and the MD trajectories were recorded every 10 ps for the following analysis. The root mean square deviation (RMSD), root mean square fluctuation (RMSF), and hydrogen bonding was measured and visually analyzed in the results section.

Estimation binding free energy via MMPBSA
Molecular mechanics/Poisson-Boltzmann surface area (MM/PBSA) [39][40][41] was applied to determine thermodynamical stability of small molecules inside the binding sites of the proteins which were computed based on the equations shown below: In formula (1), the binding free energy ( ∆G binding ) can be decomposed into five terms: internal energy ( ∆E internal ), van der Waals ( ∆E VDW ), electrostatic interaction ( ∆E elec ), and free energy of solvation which consists of polar ( ∆G GB ) and non-polar solvation free energy ( ∆G SA ). Polar solvation free energy was calculated by the GB model (igb = 2) [42] and non-polar solvation free energy was calculated based on the solvent-accessible surface area (SASA) estimated by the LCPO algorithm: ∆G SA = 0.0072 × ∆SASA [40]. All (1) the MM/GBSA free energy calculations were performed by using the MMPBSA module in the AMBER 18 package.

Statistical analysis
The results were expressed as mean ± SEM analyzed with a one-way analysis of variance (ANOVA) using Prism 8.0.1 (GraphPad, San Diego, CA, USA). P < 0.05 was considered to be significantly different.

Growth inhibition of AGS on CRC
MTT assay was used to evaluate the influence of AGS on the viability of the CRC cells. As shown in Fig. 2 (Fig. 2b), respectively, indicating that the three fractions of AGS could significantly inhibit the growth of CRC cells. These data suggested that AGS would be a good candidate for CRC prevention or therapy.

Collection for disease targets
A total of 21,572 targets related to CRC disease were obtained by searching the Genecards database (Supplementary Table S2). The top 10 "high response" genes were screened out according to the score of relevance, including BRCA2, BRCA1, TP53, MSH2, APC, MSH6, MLH1, CDH1, PTEN, and PMS2.

Prediction for candidate targets of AGS against CRC
As the Venn diagram showed in Fig. 4a, a total of 170 overlapped genes (Supplementary Table S3) were identified by matching the therapeutic target genes of CRC and target genes of AGS. The "AGS-component-target-CRC" network was built up by importing the crossover genes of AGS & CRC and the potential active components into the system (Fig. 4b).

PPI network construction and drug-disease key targets prediction
Based on the crossover targets of AGS & CRC, the PPI network was built up by the String database. As shown in Fig. 5a, there were 169 nodes and 239 edges in the network diagram and the average node degree was 2.83. The top 30 core genes were screened out, of which the node degrees of SRC, MAPK1, ESR1, HSP90AA1, and MAPK8 were greater than 12 (Fig. 5b). According to the results, the core targets were predicted to be SRC, MAPK1, ESR1, HSP90AA1, MAPK8, which had more connections than other genes.

Enrichment analysis for key targets
As a result, GO analysis showed that the numbers of BP, CC, and MF of AGS against CRC were 1079, 44, and 132, respectively (Supplementary Table S4), and the top 20 GO analyses of BP, CC, MF had been shown as graphical bubbles. According to the results, AGS would mainly participate in the biological process of steroid metabolic process (Fig. 6a) and the cytoplasm is the major reaction site (Fig. 6b) in the treatment of CRC, during which the central molecular function would probably include the steroid hormone receptor activity, nuclear receptor activity, transcription factor activity, steroid binding and endopeptidase activity (Fig. 6c). Furthermore, a total of 96 signaling pathways (Supplementary Table S5) were screened out through KEGG pathway enrichment analysis, and the top 20 signaling pathways were shown as a bar graph (Fig. 6d), among which the MAPK signaling pathway, lipid, and atherosclerosis, proteoglycans in cancer, prostate cancer, adherens junction, endocrine resistance, progesterone-mediated oocyte maturation, relaxin signaling pathway, foxO signaling pathway, and apoptosis had been proved to be the major pathways related to CRC treatment.

Molecular docking analysis
The docking analysis was performed to assess the binding effect and pattern between the active ingredients of AGS and the identified core targets. The results of SP molecular docking between SRC (PBD ID: 3G5D) and the 6 characteristic active components were shown in Table 1   Altogether, the results showed that the target proteins SRC and ESR1 had stronger docking capability with AGS than the other targets. Four compounds (11-O-galloylbergenin, 11-O-protocatechuoylbergenin, 11-O-syringylbergenin, epicatechin-3-gallate) had a lower docking score than the ligand-protein of SRC, and two compounds (bergenin, epicatechin-3-gallate) docked better than the ligand-protein of ESR1 (Table 1), which indicated that AGS had multiple ingredients and multiple targets against CRC.

Stability analysis
In this study, molecular dynamics (MD) simulation was used to analyze the stability of compounds from AGS binding multiple targets including SRC, MAPK1, ESR1, HSP90AA1, and MAPK8. We analyzed the root mean square deviation (RMSD) plot, presenting the stability of protein, and the higher of RMSD value is, the more violent the fluctuation is. On the contrary, the protein binding is stable. As shown in Fig. 8, all the systems were generated and submitted for 100 ns in MD simulations. Among all the systems, HSP90AA1-11-O-protocatechuoylbergenin complex and HSP90AA1-11-O-syringylbergenin complex (Fig. 8d) were the most stable systems, with the lowest RMSD value fluctuating stably within the range of 2 Å, which were followed by SRC-11-O-protocatechuoylbergenin complex, SRC-epicatechin-3-gallate complex (Fig. 8a) and MAPK8-11-Osyringylbergenin complex (Fig. 8e), fluctuating steadily at about 2.4 Å. The RMSD values of ESR1-bergenin complex and ESR1-epicatechin-3-gallate complex (Fig. 8c) were much higher than other systems, while these two systems had a high degree of convergence and stable fluctuations in the late simulation week, suggesting that they were also very stable in the middle and late simulation. For MAPK1-epicatechin-3-gallate complex, MAPK1-11-O-galloylbergenin complex (Fig. 8b), and MAPK8-11-O-galloylbergenin complex (Fig. 8d), the fluctuations were slightly strong during the simulation process, while the RMSD did not exceed the 0.3 Å threshold [41], indicating that there were no remarkable conformational changes.

Residue mobility analysis
To examine the structural-flexibility effect of small molecules on target proteins in the process of molecular dynamics simulation. We analyzed the root mean square fluctuation (RMSF) plot of the 10 complexes. In general, the flexibility of the target protein decreases after a compound binding to the protein, to stabilize and activate the protein [43]. As shown in Fig. 9, the RMSF values (color line) of ESR1-bergenin complex, ESR1-epicatechin-3-gallate complex, HSP90AA1-11-O-protocatechuoylbergenin complex, HSP90AA1-11-O-syringylbergenin complex, MAPK8-11-O-galloylbergenin complex, and MAPK8-11-O-syringylbergenin complex were lower than these in the apo form of corresponding proteins (black line) at multiple amino acid sequences (such as ESR1 210-250, HSP90AA1 55-60, 100-110, MAPK8 20-45, 60-70,100-115 regions). These results indicated that the fluctuations of backbone atoms in these proteins were significantly reduced after binding small molecules, and the proteins became more stable compared with apo-form protein. However, the same changes were not observed in the complexes of SRC and MAPK1, meaning that the fluctuation increased due to the binding of small molecules.

Hydrogen bonds analysis
The hydrogen bond is one of the strongest non-covalent interactions between small molecule compounds and proteins. During the molecular dynamic simulation, we analyzed the number of hydrogen bonds in all 100 ns of 10 complexes. As shown in Fig. 10, the hydrogen bonds could be seen in all complexes and the average number of hydrogen bonds were 4 in ESR1-epicatechin-3-gallate complex, HSP90AA1-11-O-protocatechuoylbergenin complex, MAPK1-11-O-galloylbergenin complex, SRC-11-O-protocatechuoylbergenin complex was 4, being much more than other complexes, which suggested that these four complexes were more stable than others. The ESR1-bergenin complex, HSP90AA1-11-O-syringylbergenin complex, MAPK1-epicatechin-3-gallate complex, MAPK8-11-O-galloylbergenin complex, MAPK8-11-O-syringylbergenin complex and SRC-epicatechin-3-gallate complex formed less than 4 hydrogen bonds in the middle and late stages of simulation, indicating that the formation of these complexes may not depend on hydrogen bonds.

Discussion
Colorectal cancer ranks third among the causes of cancer-related death and carries a huge global health burden [44]. Although the development of diagnostic and therapeutic technology has improved the survival of patients, approximately 0.8 million CRC-related deaths occur worldwide each year [3]. Chemotherapeutic agents such as ramucirumab, aflibercept, and bevacizumab have been developed to use in the clinic, while they would induce non-negligible side effects. Surgery remains the preferred treatment and 5FU is still the most commonly used drug for CRC therapy [45]. Hence, new agents with more efficacy and low side effects are urgently needed for CRC therapy.
Traditional Chinese herbal medicine (TCM) has been used for thousands of years in China, with the characteristics of definite therapeutic effects and low side effects, and has been used as an alternative therapy for cancer patients. Shi et al. found that TCM could significantly improve the disease-free survival of cancer patients, in particular with stage III patients [10], and Xu et al. found that long usage of TCM herbals could improve the survival outcomes in stages II and III CRC patients in China [9]. Meanwhile, a series of compounds derived from TCM such as baicalein (a constituent of Scutellaria baicalensis Georgi, Huangqin), curcumin (a constituent of Curcuma longa L., Jianghuang), berberine (a constituent of Coptis chinensis Franch., Huanglian) have been demonstrated to have good efficiencies in cancer treatment [46][47][48]. Therefore, TCM can be recognized as an important resource for antitumor new drug development.
In our study, we verified that the 3 parts of AGS extracts, PEAGS, EAAGS, and NBAGS showed a good anti-CRC effect against HCT-116 and SW620 cells proliferation determined by MTT assay. Literature searching indicated that the flavones (4-O-galloylbergenin, 11-O-galloylbergenin, 11-O-protocatechuoylbergenin, 11-O-syringylbergenin, bergenin, epicatechin-3-gallate, quercetin) and terpenoids (ardisiacrispin B, gallic acid, stigmasterol, stigmasterol-3-o-β-D-glucopyranoside) were the major pharmacological active ingredients of AGS [49,50]. Chang found that a high intake of flavones (such as quercetin) may reduce the risk of colon cancer [51]. Bergenin, one of the active ingredients of AGS, had been found to inhibit bladder cancer progression by activating the PPARγ/PTEN/Akt signal pathway [52]. Epicatechin-3-gallate was also found to possess a series of pharmacological and physiological properties, including induction of phase II enzymes, mediation of anti-inflammation response, regulation of cell proliferation and apoptosis effects, and prevention of tumor angiogenesis, invasion, and metastasis [53]. Numbers of studies indicated that quercetin had antitumor effects by promoting cell apoptosis, and autophagy, and inhibiting MAPK/Erk, PI3K/Akt, and NF-κB signaling pathways [54,55].
Network pharmacology, an emerging discipline that is commonly used to predict the potential targets and mechanisms based on the accumulation of evidence from big data, would provide strong support and accelerate the process of new drug development. In our study, we utilized the network pharmacological methodology to explore the pharmacological mechanism of the active components of AGS against CRC, and molecular docking and molecular dynamic simulation were used to investigate the docking pattern and capability between major active compounds of AGS and potential target proteins. The results showed that 11 active ingredients in AGS would have effects on 170 overlapping genes that acting important roles in CRC treatment. PPI network study showed that SRC, MAPK1, ESR1, and HSP90AA1 were the most correlated proteins, followed by MAPK8. SRC family kinase is a key mediator of cellular tumor-promotion genic signals linked with tumor proliferation, migration, and invasion, and SRC could be activated in CRC through various mechanisms, including the regulations of the SRC-STAT3 signaling pathway, SRC-CTNNB1, and macro-autophagy/autophagy pathways [56,57]. Therefore, SRC inhibitors have been considered ideal therapeutic agents. Mitogen-activated protein kinase (MAPK) is a kinase family that converts extracellular stimuli into various cellular responses and participates in the process of disease occurrence and development [58]. Multiple MAPK pathways are associated with the processes of mitosis, differentiation, metabolism, motility, apoptosis, and survival of eukaryotic cells [59]. MAPK1 (mitogen-activated protein kinase 1), an important member of the MAPK kinase family, includes a highly conserved serine/threonine kinase domain [60], and has been reported to be a target for microRNAs, such miR-212 [61], miR454 [62] and miR-422a [63]. MAPK8 at phosphorylation status was also found to be associated with germ cell apoptosis and redistribution of the Bcl2modifying factor [64]. ESR1, one of the estrogen receptors (ER), is activated by the sex hormone estrogen. A study showed that approximately 5% of primary tumor patients harbored ESR1 mutation which was increased to 30% ~ 40% in the metastatic cases [65]. A recent study indicated that ESR1 participated in the development and progression of CRC, leading to the inferior clinical outcome of CRC patients [66]. Hence, ESR1 was regarded to be an important therapeutic target against CRC. Heat shock proteins (HSPs) are commonly over-expressed in many kinds of tumors and are highly associated with a poor prognosis and therapy resistance [67]. HSP90 is one of the HSPs proteins, which has been found to restrain cell apoptosis through folding, stabilizing, and activating oncogenic proteins [68]. HSP90AA1, a member of the HSP90 family that is expressed extracellularly, is strongly associated with cancer cell invasion [69]. Therefore, targeting HSP90AA1 can be considered a good remedy for CRC therapy. Based on the predicted therapeutic targets and active compounds of AGS, we performed a molecular docking analysis for 6 active components and 5 collective targets to provide a rational explanation for the anti-CRC effect of AGS, and the results showed that AGS had a strong docking capability with SRC and ESR1. Furthermore, the results of molecular dynamic simulation also showed that these complexes fluctuate little during the simulation process, and the MD trajectories analysis, including RMSD, RMSF, and hydrogen bonding, showed that the compounds of AGS enhanced the stability of target proteins. GO analysis showed that the numbers of BP, CC, and MF of AGS against CRC were 1079, 44, and 132, respectively. It was concluded that AGS mainly participated in the biological process of steroid metabolic process. The major reaction site was the cytoplasm in CRC treatment, and the central molecular function mainly included steroid hormone receptor activity, nuclear receptor activity, transcription factor activity, steroid binding, and endopeptidase activity. KEGG pathway analysis showed that 96 pathways in all were included and the top 10 pathways were MAPK signaling pathway, lipid and atherosclerosis, proteoglycans in cancer, prostate cancer, adherens junction, endocrine resistance, progesterone-mediated oocyte maturation, relaxin signaling pathway, FoxO signaling pathway, apoptosis, among which the MAPK signaling pathway, lipid and atherosclerosis, proteoglycans in cancer, prostate cancer, FoxO signaling pathway, and apoptosis were commonly recognized to be cancerrelated. It is reported that the MAPK signaling pathway plays an important role in cell proliferation, and is commonly activated by its upstream growth-factor receptors, such as the epidermal growth factor which is commonly over-expressed in colorectal cancer [70]. The altered metabolism of lipids is a hallmark in many cancers. A series of lipid molecules, such as fatty acids, polar lipids, and oxylipins, can promote the development of CRC, and the lipid metabolism pathways have become the targets of CRC treatment [71]. Proteoglycans are a group of molecules that have a glycosaminoglycan chain. Previous studies showed that some proteoglycans, such as glypicans, agrin, and versican, play a key role in the development of liver cancer, and heparan sulfate proteoglycans were regarded to be the critical targets for the diagnosis and therapy of CRC [72,73].
AGS has a wide range of biological properties including anti-cancer activity. Gu found that the biotransformation product S1 from AGS had significant inhibition on 6 kinds of tumor cell lines and the potential mechanism may be related to cell cycle arrest [74]. Cyclin D1, frequently over-expressed in ESR1-mutated breast cancer [75], could activate the CDK4 and CDK6 to facilitate cell cycle progression through the G1 restriction point  [76,77]. Similarly, v-SRC (viral-SRC, a transforming protein of SRC family) was found to suppress the expression of the cyclin-dependent kinase (CDK) inhibitor p27, leading to rapid transit of the G1 phase and the expressions of CDK2, CDK4, and CDK6 [78][79][80]. He found that the components of AGS suppressed the growth of HepG2 cells by regulating the phosphorylation of ERK, JNK, and p38 in the MAPK signaling pathway [81]. The results were consistent with our findings that MAPK1 and MAPK8 were the core proteins in AGS against CRC and the MAPK signaling pathway was the one of top 10 KEGG pathways. In all, our study firstly investigates the anti-CRC of AGS, explores the potential mechanism, analyzes the docking patterns and binding capabilities between the active compounds and target proteins, as well as studies the stability of the complexes. The results suggest the potential application of AGS in colorectal cancer treatment or prevention for humans through integrating experimental evaluation and network study. The 11 potential active ingredients were showing a summary of the composition from AGS based on the existing database. Although the PEAGS showed higher cytotoxicity, the NBAGS and EAAGS also had good inhibitory potency against CRC cell growth. It was suspected that there probably were some unknown compounds with good anti-tumor activity, including the new compounds that have never been reported or excluded in the databases. Hence, it is needed to conduct a series of experiments to identify the constituent structures in PEAGS, which is our follow-up study.

Conclusion
This study showed that all the three kinds of fractions from AGS, including the n-butanol extract (NBAGS), ethyl acetate fraction (EAAGS), and petroleum ether fraction (PEAGS), could significantly inhibit the proliferation of CRC cells, with the IC 50 values of 197.24, 264.85, 15.45 µg/mL on HCT-116 cells, and 523.6, 323.59, 150.31 µg/mL on SW620 cells. Network pharmacological analysis suggested that eleven active ingredients were identified, including 4-O-galloylbergenin, 11-O-galloylbergenin, 11-O-protocatechuoylbergenin, 11-O-syringylbergenin, ardisiacrispin B, bergenin, epicatechin-3-gallate, gallic acid, quercetin, stigmasterol, stigmasterol-3-o-β-D-glucopyranoside. The PPI network showed that SRC, MAPK1, ESR1, HSP90AA1, and MAPK8 would probably be the core targets of AGS against CRC. GO analysis showed that the numbers of biological process, cellular component, and molecular function of AGS against CRC were 1079, 44, and 132, respectively, and KEGG enrichment suggested that 96 signaling pathways in all would probably be involved in AGS against CRC, among which MAPK signaling pathway, lipid and atherosclerosis, proteoglycans in cancer, prostate cancer, adherens junction would probably be the major pathways. The molecular docking study showed that the targets proteins SRC and ESR1 had strong docking capability with AGS than the other targets. The compounds of 11-O-galloylbergenin, 11-O-protocatechuoylbergenin, 11-O-syringylbergenin, epicatechin-3-gallate have higher docking scores than the ligand-protein of SRC, and the compounds of bergenin, epicatechin-3-gallate dock better than the ligand-protein of ESR1, suggesting that AGS has multiple ingredients, multiple targets, and multiple pathways against CRC. Our study for the first time investigates the anti-CRC potency of AGS, as well as uncovers the underlying mechanism. The results can probably provide valuable information for further study on the anti-CRC effect of AGS as well as its underlying mechanism.